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Abstract. We study pattern-forming instabilities in reaction-advection-diffusion systems. We de- 
velop an approach based on Lyapunov-Bloch exponents to figure out the impact of a spatially 
periodic mixing flow on the stability of a spatially homogeneous state. We deal with the flows 
periodic in space that may have arbitrary time dependence. We propose a discrete in time model, 
where reaction, advection, and diffusion act as successive operators, and show that a mixing ad- 
vection can lead to a pattern-forming instability in a two-component system where only one of the 
species is advected. Physically, this can be explained as crossing a threshold of Turing instability 
due to effective increase of one of the diffusion constants. 
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1. Introduction 



Reaction-diffusion systems is a well-established class of models describing various aspects of pat- 
tern formation far from equilibrium. Quite often pattern-forming fields are transported by fluid 
flows - examples range from a development of plankton patterns in oceanic flow [3j to chemical 
reactions in microchannels Q. By incorporating the flow in the model, one arrives at reaction- 
advection-diffusion models, with a much richer variety of possible phenomena (see, e.g., a recent 
review paper lfT3~l0 . In Ref. [XO it was demonstrated that a differential flow (where only some 
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components are advected) may lead to new instabilities in the system, in particular to new convec- 
tive instabilities ifToll . This study was restricted to a planar geometry, but later in [j4l [2j such an 
instability was demonstrated for a circular geometry (Couette flow) as well. As has been recently 
shown, a simple shear flow is able to destabilize the spatially homogeneous state lfT5l . A reason- 
able question to ask is whether a similar destabilization effect can be found for mixing flows. One 
might intuitively expect that as mixing smears spatial nonuniformities, it results in stabilization 
of a spatially homogeneous state (like it happens if the reaction is chaotic |[T2l ). In this paper 
we particularly address this question by presenting an approach to study pattern-forming instabili- 
ties in periodic in space mixing flows. It is based on the calculation of Lyapunov-Bloch exponents 
(cf. 0), and provides an efficient tool for finding mostly unstable patterns. As a particular example 
we consider the effect of a mixing advection on a general two-dimensional reaction-diffusion sys- 
tem capable of Turing instability and demonstrate that instability can be induced by an advection 
of one component of the reaction. 

2. Model formulation 
2.1. Continuous-time model 

A variety of situations in biological and chemical contexts can be described by the dynamics of 
two interacting species - an activator and an inhibitor. A popular model, e.g., is the Brussela- 
tor [|6l or its modifications. For spatially distributed fields, the dynamics also includes diffusion 
terms (molecular diffusion for chemical systems or irregular mobility in biological applications) 
and advection due to an imposed macroscopic velocity field V(r, t). We assume the latter to be 
incompressible, furthermore we normalize time by the characteristic advection time. In this paper 
we are interested in stability properties of a steady homogeneous distribution. Denoting the devi- 
ations of the concentrations of the species from this steady state with P, Q, we arrive at a linear 
reaction-advection-diffusion system governed by a couple of dimensionless equations 



Here P and Q are deviations of chemical concentrations from the steady state, Dp and Dq are their 
corresponding diffusivities; a, b, c, d are parameters of the kinetics. We assume that concentrations 
do not influence the flow, so that the velocity field V(r, t) does not depend on P and Q. 

In the absence of advection the problem reduces to the classical reaction-diffusion model (see, 
e.g., Ref. flU), where two principal instabilities are a spatially homogeneous Hopf bifurcation and 
a Turing instability |fl4l|. In this paper our main interest is in the case where Turing instability is 
dominant, what requires in particular that the diffusion constants Dp and Dq are different. Our 
goal is to describe unstable modes in the presence of advection term in (12.11) . (12.21) . This can be 
done numerically, and for an effective calculation we formulate a discrete in time model of the 
reaction-advection-diffusion system above. 




D P V 2 P + aP + bQ, 



D Q V 2 Q + cP + dQ. 



(2.1) 



(2.2) 
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2.2. Discrete in time model 

A typical experimental realization of a two-dimensional mixing flow V(r, t) is a flow periodic 
in space and time, as realized e.g. in experiments IfTOll . Qualitative behavior of the system can 
be understood from a simple problem, where advection is modeled by a two-dimensional time- 
dependent flow which is 27r-periodic in space and T-periodic in time and first introduced by An- 
tonsen et al. HI 



V = e x U x f(t) sm[y + x (t)} + e y U y [l - f(t)] sm[x + 6 y (t)}, 



(2.3) 



The function f(t) describes switching between two shear flows in x and y directions within time 
intervals T\ and T 2 with the amplitudes U x , U y , respectively 



/(*) 



1, > t > T u 
0, Ti < t < T 2 . 



(2.4) 



The advantage of such a setup is that the trajectories of particles in the flow (12.31 ), (12.41 ) can be 
obtained explicitly. The transformation of particle coordinates due to the advection during one 
time interval T = T\ + T 2 is governed by a map 



x = x + U x Ti sm(y + 6 X ) , 
y = y + U y T 2 sm(x + 9 y ) , 



(2.5) 
(2.6) 



where x = x(t n ), y = y(t n ), x = x(t n + T), y = y(t n + T), and n is the number of iteration. 
The phase space dynamics of the area-preserving map (12.51) , (12.61 ), restricted to the basic domain of 
periodicity < x, y < 2n is demonstrated in Fig. [TJ It is typical of Hamiltonian flows with regular 
islands and domains of chaotic behavior. With growth of advection rate (U x , U y ), domains of the 
quasiperiodic dynamics are gradually superseded by the regions of the chaotic dynamics. 





Figure 1: Phase portraits of the map (12.51) . (12.61) for periodic driving, Ti = T 2 = 0.5, 9 X = 9 y — 0. 
The parameters U x = U y — 1.5 (left) and U x — U y = 3.0 (right). 
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While considering of the reaction-advection-diffusion system (12.11 ), (12.21 ) with the flow ( 12.31) , 
(12.41 ) remains a complex computational problem, there is a possibility to simplify the analysis, 
reducing the system to a discrete-in-time model. In the vein of the approach used in Ref. (8]| we 
consider a model situation where the action of advection, diffusion, and reaction is separated in 
time. We assume, that within each time interval T both scalar fields P and Q evolve in three stages, 
corresponding to advection, diffusion, and reaction. 

During the first stage, both scalars evolve independent of each other according to pure advec- 
tion. For a scalar density, this process can be described in terms of the Frobenius-Perron operator 
La for the map (12.51) . (12.61) . Then, on the next step, the diffusion operator = exp(TL>jV 2 ) is 
applied. Finally, the fields are subject to reaction according to the system 

dP 

— = aP + bQ, (2.7) 
at 

^ = cP + dQ. (2.8) 
at 

Evolution during time T of this system leads to a linear transformation 

p\ /p\ / (A--a)e A + r -(A + -a)e A - r b{e x - T -e x + T ) 

Q j = LR I, Q ) = {X--a){\+-a)&+ T -e>- T ) (a-A + ~a)e x - T 1 ( ^ 1 ^ 

7 V 7 \ 6(A_-A+) A_-A+ 




where the exponents 



a + d ,\/{a - d) 2 + Abe 
A+ = — - — ± 



2 2 

are assumed to be real, in accordance with our choice of absence of Hopf bifurcation. The reaction- 
advection-diffusion propagator over one time interval is given by the product L r L Di La of the 
operators. The goal of the stability analysis is to find unstable eigenvalues and eigenmodes of this 
operator. 



3. Bloch ansatz 

Although the basic flow V(r, t) is periodic in space, a perturbation of the field has not to be peri- 
odic. A general perturbation should be taken in the Bloch form 

P{x,y,t n ) = e iK * x+iK »v${x,y,t n ), Q(x,y,t n ) = e iK * x+iK ^{x,y,Q, (3.1) 

where the functions $, ^ are 2n -periodic in space and additional parameters k x , n y stand for 
quasimomenta. Since the exponential factor e lK * x + tK yy has a period of unit with respect to k x , 
K y , we consider a symmetric interval of independent values n x , n y £ [— |, |]. Then, because of 
periodicity of functions $, ^, the solutions can be represented as Fourier series 

*(x,y,t n ) = ^Utn)e i{lx+my \ V(x,y,t n ) = J^^iQ e iilx+my) . (3.2) 

l,m l,m 
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Of all operators L A , L D , and L R the most nontrivial is the advection operator L A . We derive 
its Fourier representation in appendix [A] The resulting transformation of the Fourier components 
during the advection and diffusion stages reads 

e -[( i+Kl ) 2 + K%) 2 ]^^ fc ^ M) (33) 

p,g 

e -[(l + . x ? + (m +Ky r ]DQ T J- g lmpq ^ (3.4) 

p,q 

e ^<i-™)o. e -i<?-Wvj q _ m [ Ux T 1 (p + K x )]J p ^[U y T 2 (m + K y )}, 

where J m (z) is the Bessel function of the first kind. These components further interact according 
to the reaction stage (12.91) (where one should replace P — >• <p, Q — ► tp). The resulting model is a 
composition of (1331) . (l34l) . and (f2T9l) . 

4. Results of the stability analysis 

In this section we apply the proposed model to study the influence of the advection on pattern 
formation in a reaction-advection-diffusion system. To some extent, the influence of advection 
in an advection-diffusion system can be understood from the idea of effective diffusion: mix- 
ing effectively increases the diffusion constant. Therefore one can expect that the dynamics of a 
reaction-advection-diffusion system is similar to that of an reaction-diffusion system with larger 
diffusion constants. The mostly interesting point is that in the system under consideration there are 
two coupled species, and the Turing instability is caused by a difference in diffusivities of species. 
Although advection can effectively change the diffusion constants of the species, it is not clear how 
this difference will be influenced by advection. However, the situation becomes much more trans- 
parent if only one specie is advected. Then an advection, contributing to its effective diffusion, 
may increase or decrease the difference of diffusion constants, thus enhancing or suppressing the 
instability. 

Below we focus on a situation when only one species, namely that of higher diffusion constant, 
is advected. We set the parameters b = 8, d = —9, D P = 0.0025, Dq = 0.0075. So, we assume 
that the mobility of the "activator" P is relatively low, and it is not advected at all. 

To perform a linear stability analysis of spatially homogeneous states with respect to inhomo- 
geneous perturbations, we apply the method of (transversal) Lyapunov exponents (LE). We use the 
usual method for estimation of the largest LE of mappings (see, e.g., Ref. (HI). We start with an 
arbitrary initial distributions for <p, ip, with vanishing spatial average, and iterate the mapping, per- 
forming a renormalization of the linear fields. The averaging of the logarithm of the normalization 
factors yields the LE. Note that this method can be equally well applied to both regular and irreg- 
ular flows (12.31) , (12.41) . Here, however, we focus on the simplest case of time-periodic advection, 
when T 1 = T 2 = 0.5, T = 1, 9 X = 0,6 y = 0, and also put U x = U y = U. 

We choose the parameters of the reaction in such a way that the homogeneous solution is stable 
in the absence of advection, and then switch on mixing of the specie Q. The dependence of the 
largest LE on the advection rate is presented in Fig. [21 One can see that the impact of advection 



LDpLA^lm 
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Figure 2: Largest Lyapunov exponent as a function of advection rate U. The parameters are 

a = 5.35, c = —6.35, k x = 0.5, k v = 0.5. 



results in the growth of the largest LE, which becomes positive at U cr ~ 0.18. So, this example 
clearly demonstrates that mixing can play a destabilizing role. 

Remarkably, the quasimomenta in the Bloch ansatz (13.11) are essential in the stability analysis. 
Below we present three examples where mostly unstable modes correspond to different values of 
quasimomenta. 

We start with the parameters of Fig. [2l Near the threshold value of advection rate we present 
the dependence of the largest LE on the quasimomenta k x , k v , see Fig. [3] (left panel). This LE 
reaches its maximum value at k x = ±0.5, K y = ±0.5, which indicates that the unstable patterns 
have a "chess-board" structure with respect to the periodicity of the original flow. A typical pattern 
of the field $ is presented in Fig.[3](right panel). 
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Figure 3: Largest Lyapunov exponent as a function of 
at a point of its maximum (right) evaluated at a = 5.35, c 
correspond to k x = ±0.5, K y = ±0.5. 



y (left) and a corresponding pattern 
= -6.35, U = 0.18. Maxima of LE 



Another set of parameters is presented in Fig. @J Here, the maximum of LE corresponds to 
k x = 0, Ky = and the periodicity of the pattern is the same as the periodicity of the imposed flow. 
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Figure 4: Largest Lyapunov exponent as a function of k x , k v (left) and a corresponding pattern at a 
point of its maximum (right) evaluated at a = 3.7, c = —4.7, U = 0.77. Maxima of LE correspond 

tO K x = 0, K y = 0. 

Finally, in Fig. [5] we show the case where the maximum of the largest LE corresponds to the 
points k x = 0, K y = ±0.4 and k x = ±0.4. This is an example of a nontrivial situation where the 
periodicities of the patterns in x- and y- directions are not the same and are not identical to the 
periodicity of the flow. We stress that in all figures above we have shown a linear mode with the 
largest growth rate. On a nonlinear stage (which is beyond the scope of this paper) a pattern of the 
flow may significantly deviate from the linear one. 
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Figure 5: Largest Lyapunov exponent as a function of k x , k v (left) and a corresponding pattern at a 
point of its maximum (right) evaluated at a — 3.1, c = —4.1, £7 = 3.5. Maxima of LE correspond 
to k x = 0, k v = ±0.4 and k x = ±0.4, k v = 0. 



5. Conclusion 

The main goal of this paper is to develop an approach, based on Lyapunov-Bloch exponents, for 
an analysis of pattern-forming instabilities in reaction-advection-diffusion systems. It is applicable 
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to periodic in space flows that may have arbitrary time dependence. We demonstrated the method 
using a discrete in time model, where reaction, advection, and diffusion act as successive operators. 
This enormously simplifies the calculations, while yielding a qualitatively correct picture of the 
process. For an exact quantitative analysis of the full system one has to apply the Lyapunov-Bloch 
ansatz to the full equations. 

We have demonstrated that a mixing advection of one of the species may lead to a pattern- 
forming instability. Physically, this can be explained as crossing a threshold of Turing instability 
due to effective increase of one of the diffusion constants. Of course, mixing can play also a 
stabilizing role, suppressing spatially inhomogeneous perturbations (see, e.g., Ref. lfT2l for such 
an analysis of stabilizing role of advection for chaotic in time reaction). Nonlinear patterns beyond 
the transition deserve further investigation, which goes beyond the scope of this paper. 
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A Advection-diffusion map with the Bloch ansatz 

Consider the transformation of a scalar field <fi due to advection, according to the map (12.51) . (12.61) . 
We assume that 

P(x,y,t n ) = e^ x+iK y y $(x,y,t n ) 

and 

l,m 

As the map is required to be area-preserving, we can write for any iteration P dxdy = P dxdy 
(where P = P(x, y), x = x(t n + T),y = y{t n + T)) or 

e iK x x+i Ky y <| ^dy = e iKxX+™yV $ dxdy. (Al) 

According to (12.51 ), during the first time interval 2\ we have x = x + U x Ti sin(y + 9 X ), y — y, 
and therefore (IA1I) gives us 

e iK x U x T lS in(y+d x ) q d - d - = $ dxdy 
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Taking this result into account, for 0; m = 4>i m {t n +T X ) we successively obtain 

4> lm = JJ § e - i ( te +™0 dxdy = j^- Jj \ e -™*U x T lS in(y+O x ) g-iflx+mj?) d%dy 

— \ J J e -iK x U x T! sm(y+0 x ) e -i[lx+lU x Ti sm(y+6 x )+my] dxd y 
= I f e -iKxU x TxBin(y+O x ) e -i[lx+lU x Ti sin(y+0 x )+my] e i(px+qy) dxd y 

^ 7r ^ p,q 

2tt 2tt 

= _L_^0 pg [ e i(p - l)x dx f e i{q - m)y e- iU * Tl{l+K * )sin( - y+8x) dy 

71 P-9 

= ^vY.^hh = -^^Cp 9 27r^2 7 re-^-^J 9 _ m [^T 1 (m + ^)] 

= ^e-^- m ^J,_ m [f/ :c T 1 (m + K x .)]^, (A2) 
<? 

in virtue of the integrals 

Mv—l)x j_ o_ 



h = J e i{p - l)x dx = 27r5 p 
o 

h = J e i to- m h e -™* n ^ ) ^W* ) dy = e- i ( q - m)e *2n J q _ m [U x Tx{l + k x )}, 
n 

where the Bessel function of the first kind appears J m (z) = ^ J Q 27r e im< ' e -;i2;sm ^<iC. 

During the second time interval T 2 we have x = x, y = y + U y T 2 sin(x + #,,), which according 
to (IA1 1) leads to an equality 

e i KH c/ H r 2 sin(x+e H ) ^ = $ dxdy _ 
From an analogous consideration we obtain for (f> lm = Jm = 4>i m (t n + Ti) = 4>i m (t n + Ti + T 2 ) = 

<&m = e-^H-i[*W(m + S )]0 pm . (A3) 



Thus, combining (|A2I) . (I A3 1) we obtain the map describing action of advection operator La 
within time interval T 

L A <f>im = Yl e-^- m)e * e-^H-^Tifr + «x)] Jp-i[^ 2 (m + ^)]0 pg . (A4) 

P.9 
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Now consider the action of the diffusion operator L D = exp(7\DV 2 ), describing diffusion 
spreading of the passive scalar with the diffusivity D within time interval T. Because this operator 
is diagonal, it acts on each Fourier component independently: 

Wim = e-^ + ^ m+ ^ DT <p lm . (A5) 

The action of both advection (IA4I) and diffusion (IA5I) operators leads to the following transfor- 
mation of the Fourier amplitudes 

L D L A cj> lm = e-^ 2+ ^ +K ^ DT J29im Pg <f> Pg , (A6) 

p,g 

9i mpg = e- i b- m '> *e- i WyJ q - Tn [U x T 1 (p + K x )}J p ^[U y T 2 (m + n y )}. 
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